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Abstract. We compare the general effective theory of one-body dark matter nucleon interac¬ 
tions to current direct detection experiments in a global multidimensional statistical analysis. 
We derive exclusion limits on the 28 isoscalar and isovector coupling constants of the theory, 
and show that current data place interesting constraints on dark matter-nucleon interaction 
operators usually neglected in this context. We characterize the interference patterns that can 
arise in dark matter direct detection from pairs of dark matter-nucleon interaction operators, 
or from isoscalar and isovector components of the same operator. We find that commonly 
neglected destructive interference effects weaken standard direct detection exclusion limits 
by up to one order of magnitude in the coupling constants. 
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1 Introduction 

The simultaneous operation of dark matter direct detection experiments exploiting comple¬ 
mentary techniques, and target materials, raise the issue of how to effectively interpret the 
information gathered on dark matter [1-9]. 

The currently favored strategy in the analysis of dark matter direct detection experi¬ 
ments relies on drastic simplifying assumptions regarding the underlying dark matter-nucleon 
interaction [10, 11]. For instance, the dark matter-nucleon interaction is often assumed to be 
independent of the momentum transferred in the scattering, and of the dark matter-nucleus 
relative velocity. Moreover, the possible interference of different dark matter-nucleon inter¬ 
actions is commonly neglected in the data analysis. At the same time, distinct experiments 
are usually interpreted separately, even when they are compatible and, most importantly, 
complementary. 

Though this approach is a well motivated first approximation, the substantial progress 
recently made in the field, and the efforts planned for the next decade, motivate the explo¬ 
ration of more sophisticated strategies. 
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In this work we explore two extensions of the standard approach to dark matter direct 
detection. Firstly, we interpret current dark matter direct detection experiments without 
assuming the knowledge of the dark matter-nucleon interaction a priori. In other words, 
we do not artificially restrict the analysis to constant spin-independent or spin-dependent 
interactions. In contrast, we interpret current observations within the general effective the¬ 
ory of one-body dark matter nucleon interactions, which predicts 28 independent isoscalar 
and isovector dark matter-nucleon interaction types. Secondly, we interpret current direct 
detection experiments globally [12, 13], i.e. varying all correlated coupling constants simul¬ 
taneously, and considering several direct detection experiments in a single multidimensional 
Likelihood analysis. 

Our approach allows to characterize the interference patterns arising in dark matter 
direct detection from pairs of dark matter-nucleon interaction operators, or from isoscalar 
and isovector components of a given operator. We refer to these interference patterns as 
multi-interaction interference effects. A key result of this work is to quantify the impact 
of multi-interaction interference effects on the calculation of dark matter direct detection 
exclusion limits. 

Effective theories for dark matter-nucleon interactions were already explored in the 
context of dark matter direct detection in [13-33], in the analysis of neutrino telescope 
observations in [34-38] , and in the context of helioseismology in [39-42] . The prospects for 
dark matter direct detection in general effective theories were first investigated in [26-28], 
and later in [32] . Ref. [26] also quantifies the bias in the dark matter particle mass and cross- 
section reconstruction induced by fitting observations assuming constant dark matter-nucleon 
interactions, when dark matter interacts with nucleons differently. The importance of multi¬ 
interaction interference effects in dark matter direct detection was first noticed in [13, 43, 44], 
and recently also in [32]. 

The paper is organized as follows. In Sec. 2 we review the non-relativistic effective 
theory of dark matter scattering from nucleons and nuclei. Sec. 3 is devoted to the statistical 
framework adopted in the analysis of the dark matter direct detection experiments introduced 
in Sec. 4. We present our results in Sec. 5, and conclude in Sec. 6. Important equations are 
listed in Appendix A. 

2 Dark matter-nucleus scattering in effective theories 

In this section we review the non-relativistic effective theory of dark matter scattering from 
nucleons and nuclei. 

2.1 Interaction Hamiltonian 

Assuming one-body dark matter-nucleon interactions mediated by a heavy spin-1 or spin-0 
particle, the most general Hamiltonian density for dark matter-nucleus scattering is given 

by [17] 

= E E (2-1) 

i=l T=0,1 k 

It is the sum of A terms, one for each nucleon in the target nucleus. The A terms in the sum 
are linear combinations of 14 Galilean invariant quantum mechanical interaction operators 

fc = 1, 3,..., 15. They are listed in Tab. 1. 
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" (i) 

The 14 operators Oj^' are constructed from the momentum transfer operator q, the 
relative transverse velocity operator v-*-, and the dark matter particle and nucleon spin 
operators, and Stv, respectively [17]. In Eq. (2.1) = 1 is the identity in isospin space, 

= Ts is the third Pauli matrix, and mjsf is the nucleon mass. The isoscalar and isovector 
coupling constants, respectively and c^, are related to the coupling constants for protons 
and neutrons by c| = (c° + c^)/2, and = (c° — c\)/2. They have dimension mass to the 
power of —2. 

The Hamiltonian density in Eq. (2.1) admits the following coordinate space representa¬ 
tion [17] 


'WT(r) = 


{ A A 

^ % <5(r - Fi) + • Ui 6{r - r^) 

i=l i=l 

A ^ r 

+ • 7^- i^r5(r - ri) - *(^(r - ri)'^r 

^ 2m]v 

2 = 1 

A ^ 

+ X ai5(r - Fi) + (5(r - Fj) (Tj X 

IrriN 


2 = 1 


A^) 


( 2 . 2 ) 


where ai denotes the set of three Pauli matrices representing the spin operator of the ith- 
nucleon in the target nucleus. Fj is the rth-nucleon position in the nucleus center of mass 
frame, and 


% — c\ + i( X Cg + Cg -)- c[i 
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F — - 
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qr,cr q^cr qq^r 

— Cg + zSx Ci2 - — X Sx Ci3 - I — — ■ Sx Ci5 
777tV 777 tv 777 TV 777tv 


(2.3) 


In the last expression v;^ = v-*- — where is an operator acting on the ith-nucleon 
space coordinate. Explicit coordinate space representations for the operators q, and 
can be found in [37]. 

Eq. (2.2) shows that dark matter couples to the constituent nucleons through the nuclear 
vector charge and spin current (first line), the nuclear convection current (second line), and 
the nuclear spin-velocity current (last line). In Eq. (2.2) we have omitted the nuclear axial 
charge, since it does not contribute to dark matter-nucleus scattering cross-sections when 
nuclear ground states are eigenstates of P and CP, as it is commonly assumed for target 
nuclei in dark matter direct detection. 
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^3 = iSN ■ X V^) 

O 4 = • Sjv 

^5 = iS^ ■ X v^) 

= (Sx • ^) (Sjv • ^) 
dr = Sn-v^ 
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a9 = iS^- S^x-^ 
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Table 1. Non-relativistic quantum mechanical operators defining the general effective theory of one- 
body dark matter-nucleon interactions. Because of the nucleon mass, tojv, in the equations above all 
operators have the same mass dimension. For simplicity, here and in the next sections we omit the 
nucleon index (i) adopted in Sec. 2. 


2.2 Transition amplitude 

We use Eq. (2.2) to calculate (|Ad at/jP) spins, i-e. the square modulus of the amplitude for 
dark matter-nucleus scattering. Assuming that nuclear ground states are eigenstates of P 
and CP, and averaging (summing) over initial (final) spin configurations, one finds^ 
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dvr 


NRj /spins 


2J + 


lE 
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'N 
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Rl 




m 


N 


(y) 


(2.4) 


The index k extends over the nuclear response functions W^'^'(y) defined below in Eq. (2.7). 
Assuming the harmonic oscillator basis for single-nucleon states, y = ( 6 g/ 2 )^ with 

b = y^ 41 . 467 /( 45 A-V 3 _ 25 A- 2 / 3 ) fm. ( 2 . 5 ) 

The 8 dark matter response functions K^' in Eq. (2.4) depend on matrix elements of the 
operators in (2.3), and are listed in Appendix A. They are functions of the dark matter 
particle spin, of /m?^, and of 

Vt^ = v‘^-q'^/dyT): (2-6) 


where yx and v are the reduced dark matter-nucleus mass and the dark matter-nucleus 
relative velocity, respectively. The 8 nuclear response functions (y) in Eq. (2.4) are 
quadratic in nuclear matrix elements, and defined as follows 

W^AB(y) = ^priq) \\J,T,Mt){J,T,Mt\\ BL-r'{q) \\J,T,Mt), 

L 

^Strictly speaking Mnr is not an amplitude, as with our definition of cj it has dimension mass 
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(2.7) 











where | J, T, Mt) represents a nuclear state of spin J, isospin T, and isospin magnetic quantum 
number Mt- The reduction operation, i.e. || • ||, is done via the Wigner-Eckart theorem. In 
Eq. (2.7), A and B can each be one of the following nuclear response operators 
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( 2 . 8 ) 


where Mf^(g'ri) = jL(gri)Y^^(nrJ , and MiMiq^i) = jLiqri)YLMAA- The vector spher¬ 
ical harmonics, Y^j^(flrj), are defined in terms of Clebsch-Gordan coefficients and scalar 
spherical harmonics: 


= Y.{L'mlX\L'lLM)YL^^{nr^)ex , (2.9) 

mX 


where is a spherical unit vector basis. In Eq. (2.4), we use the notation W^' (y) = W^Av) 
for A = B. 

The 6 nuclear response operators in (2.8) arise from a multipole expansion of the nuclear 
charge and currents in Eq. (2.2). The multipole index L must be less than 2J. In our analysis 
of current direct detection experiments, we adopt the nuclear response functions of Ref. [21]. 
Nuclear response functions for dark matter capture by the Sun can be found in [37] in analytic 
form. 


2.3 Scattering rate 

In a dark matter direct detection experiment, the expected differential rate of scattering 
events per unit time and per unit detector mass is given by 


d7^ 

dEji 
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dB-T 




f{v + Veit)) 
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i(<?) 
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( 2 . 10 ) 


where is the dark matter particle mass, ^t is the mass fraction of the nucleus T in the 
target material, and is the dark matter density in the solar neighborhood. In Eq. (2.10) 
the integral denotes an average over the local dark matter velocity distribution, /, in the 
galactic rest frame boosted to the detector frame, nmin(<?) = qApT is the minimum velocity 
required to transfer a momentum q from the target nucleus to the dark matter particle, and 
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Parameter 

Type 

Prior range 

Prior type 

Reference 

logio (clm^) 

model parameter 

[-5,1] 

log-prior 

- 

Iogio(c5m2) 

model parameter 

[- 2 , 6 ] 

log-prior 

- 

logio(cjm 2 ) 

model parameter 

[-3,3] 

log-prior 

- 

logio(c^m 2 ) 

model parameter 

[- 2 , 6 ] 

log-prior 

- 

logio(cgm 2 ) 

model parameter 

[- 2 , 6 ] 

log-prior 

- 

logio(c^m 2 ) 

model parameter 

[- 2 , 6 ] 

log-prior 

- 

logio(4m2) 

model parameter 

[- 2 , 6 ] 

log-prior 

- 

logioCcgm^) 

model parameter 

[- 2 , 6 ] 

log-prior 

- 

logio(cio"i?) 

model parameter 

[- 2 , 6 ] 

log-prior 

- 

Iogio(clim 2 ) 

model parameter 

[-4,4] 

log-prior 

- 

logio(ci 2 "iS) 

model parameter 

[- 2 , 6 ] 

log-prior 

- 

loglo(cl3"i?) 

model parameter 

[- 2 , 6 ] 

log-prior 

- 

Iogio(ci4rn?) 

model parameter 

[- 2 , 6 ] 

log-prior 

- 

Iogio(cl5m2) 

model parameter 

[-1,7] 

log-prior 

- 

logio("ix/GeV) 

model parameter 

[0.5,4] 

log-prior 

- 

^Xe 

nuisance 

[0.78,0.86] 

Gaussian, a = 0.04 

[45] 

OCOUPP 

nuisance 

[0.13,0.17] 

Gaussian, a = 0.02 

[ 8 ] 

OPICASSO 

nuisance 

[2.5, 7.5] 

Gaussian, a = 2.50 

[9] 


Table 2. List of model and nuisance parameters with corresponding prior types and ranges. The 
nuisance parameters ^xe, acoupp; and opicASSO are introduced in Sec. 4 to model various types 
of detector uncertainties. Following [21], we have expressed the coupling constants in units of 
ml = (246.2 GeV)^. 


Ve{t) is the time-dependent Earth velocity in the galactic rest frame. Here we consider a 
Maxwell-Boltzmann distribution f{v -|- Ve{t)) oc exp(—IE-|- Ee(t)p/uo) truncated at the local 
escape velocity Vesc = 554 km s“^, and with vq = 220 km s“^. 


3 Statistical framework 


We compare the effective theory of dark matter-nucleon isoscalar and isovector interactions 
reviewed in Sec. 2 to current dark matter direct detection experiments in a global multidi¬ 
mensional statistical analysis. If not otherwise specified, for a given dataset d we assume the 
Likelihood function 


- - k - kin 


+ fiB 

k 


(3.1) 


In Eq. (3.1) /i 5 (m^,c, 77 ) is the number of predicted scattering events at a given point in 
parameter space, k is the number of observed recoils in a given experiment, and /is is the 
corresponding number of expected background events. The arrays c and rj represent the 28 
coupling constants of the theory and the nuisance parameters in Tab. 2, respectively. When 
the error us on the number of expected background events /is is known, we analytically 
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marginalize Eq. (3.1) over fiB obtaining in the limit as ^ 


- ln£eff(d|mx, c, 77 ) ~ c, r?) + (2 - k) In usim^, c, rf) + 


— In ' 


(/c^ - k) 2 


(Tb + 


- -^ctb 


+ constant. 


(3.2) 


The constant term in the last line of Eq. (3.2) is an arbitrary normalization, that we choose 
such that InTeff = 0 iov fj,s + kB = k. 

We construct approximate 2D frequentist confidence intervals for pairs of model param¬ 
eters from an effective chi-square defined as = —2 In £prnf/T mav , where £max is the 

maximum Likelihood, and Tprof is the 2D profile Likelihood: 

>Cprof(d|6»i,6'2) oc max £efr(d|©). ( 3 . 3 ) 

83,■■■, 9 m 

In Eq. (3.3), Q = {01,02, ■■■, 0m) = ("^x; c, 77 ). Wilks’ theorem guarantees that under certain 
regularity conditions the distribution of AXgfj converges to a chi-square distribution with 2 
degrees of freedom [46] . 

Here we mainly restrict the analysis to a frequentist approach in order to avoid the 
volume effects found in [13] studying a smaller sample of dark matter-nucleon interaction 
operators. However, in computing exclusion limits that do not involve marginalization, we 
also adopt a Bayesian approach, since it requires a smaller number of Likelihood evaluations. 
Within the Bayesian approach, we express exclusion limits on the coupling constants in 
terms of x% credible regions. Credible regions are defined as portions of the parameter space 
containing x% of the total posterior probability, and such that the posterior probability 
density function ^(0171) oc £efr(d|0)vr(©) at any point 0 inside the credible region is larger 
than at any external point. 

The function 7 r( 0 ) is the prior probability density function. We assume a uniform prior 
probability density function for the logarithm of and cl, with r = 0,1 and 7 = 1,3,...,15. 
This choice allows us to sample the Likelihood function covering several orders of magnitude 
in all directions in parameter space. Eor the nuisance parameters we assume Gaussian priors, 
as shown in Tab. 2. 

We sample the multidimensional Likelihood function (3.2) using the Multinest pro¬ 
gram [47-49]. We run Multi nest with parameters set to nuve = 20000 and tol = 10“^. 
Einally, we use our own routines to evalnate scattering rates and the Likelihood function, 
Getplots [50] to calculate the profile Likelihood, and Matlab to produce the figures. 


4 Datasets 

Here we list the experiments considered in the analysis. Eor each experiment we provide 
information needed to evaluate the corresponding Likelihood function. We refer to [13] 
for further details. In the analysis, we convolve Eq. (2.10) with a Gaussian probability 
distribution function of standard deviation cj to account for finite energy resolution effects. 

4.1 CDMS-Ge 

We consider data from the final exposure of the GDMS H experiment, corresponding to an 
exposure of 612 kg-days [1]. In the data analysis the GDMS collaboration has observed 
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k = 2 candidate events within the signal region 10 - 100 keV, with an expected number of 
background events in the same energy interval equal to = 0.9 ± 0.3. Here we assume a 
maximum experimental efficiency of 32% at 20 keV, linearly decreasing to 20% at 10 keV and 
at 100 keV. In addition, we consider a Gaussian energy resolution with standard deviation 
a = 0 . 2 . 

4.2 CDMS low threshold 

In the 2 keV threshold analysis performed by the CDMS collaboration, the T1Z5 germanium 
detector has observed k = 36 events in the 2-20 keV region [2]. Zero-charge, surface, 
and bulk events, can explain 75% of the observed counts [2]. In computing the Likelihood 
function, we therefore consider hb = 36 x 0.75, with cts = 0. We assume the efficiency 
reported in the inset of Fig. 1 in [2], a Gaussian energy resolution with energy dependent 
standard deviation a = y^O. 2932 -|- (0.05675^/?)^, and an exposure for T1Z5 of 241/8 kg-days. 

4.3 SuperCDMS 

The SuperCDMS experiment has collected data using 7 germanium detectors, with a total 
exposure of 577 kg-days [3]. Tab. 1 in [3] provides details regarding the number of recorded 
events, as well as the number of expected background events for each detector separately. In 
our analysis we consider data from 5 detectors only, hence neglecting T5Z2 and T5Z3, since 
for these detectors the number of observed counts is significantly larger than the number of 
expected background events. In the calculations we consider a Gaussian energy resolution 
with a = 0.3, and the detector efficiency of Fig. 1 in Ref. [3]. We assume an average exposure 
per detector of 577/7 kg-days. In the Likelihood function we set ^b and cjb as in Tab. 1 of 
Ref. [3]. 

4.4 CDMSlite 

The SuperCDMS experiment has also collected data in the low ionization threshold operating 
mode (i.e. CDMSlite) during 10 live days of dark matter search [4]. The nuclear recoil energy 
threshold for this experimental configuration is 170 eVee, corresponding to 841 eVnr- The 
average rate of nuclear recoils in the CDMSlite detector is 5.2 ± 1 counts/keVee/kg-day 
between 0.2 and 1 keVee, and 2.9 ± 0.3 counts/keVee/kg-day between 2 and 7 keVee. For 
these data we assume the Likelihood function 

- InTcDMSlite = 2 ®h('^[ 0 . 2 , 1 ] “ 5.2)(7?.[o,2,l] - 5.2)^ 

+ ^0H(7^[2,7] - 2.9)(7^[2,7] - 2.9)2/0.32 (4.1) 

where 7^[o.2,i] and 7^[2,7] represent the average rates between 0.2 and 1 keVee, and between 
2 and 7 keVee, respectively. We calculate the expected average count rate in the CDMSlite 
detector using the efficiency reported in the inset of Fig. 1 in [4], and assuming perfect energy 
resolution. 

4.5 XENONIOO 

We use the data presented by the XENONIOO collaboration in [5], corresponding to an 
effective exposure of 34x224.6 kg-days. XENONIOO has observed k = 2 candidate signal 
events in the pre-defined nuclear recoil energy range 6.6 - 30.5 keVnr, corresponding to 3 - 
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30 photoelectrons. The number of expected background events in the same interval is 1.0 ± 
0.2. The differential spectrum of expected photoelectrons SI in XENONIOO is given by 


d7^ 


+ 00 ^oo 

T(5'l) Gauss(51|n, -vAto'pmt) / dER Poiss(n|z^(E/j)) 

n=l Jo 


dn 

dEji 


(4.2) 


In Eq. (4.2), n is the number of photoelectrons actually produced from a recoil energy Er. 
Eor n we assume a Poisson distribution of mean ^{Er) = EiiLef[{E]i)LySnr/See, where 
Ly = 2.28 ± 0.04 PE/keVee is the light yield, See = 0.58 and ^nr = 0.95 are the electric field 
scintillation quenching factors for electron and nuclear recoils, and Lef[{E^) is the energy 
dependent scintillation efficiency. We assume for LesiEn) the following form 

( Les{Eji) for Eji/keYnr > 3 

LesiEn) = < 

[ max{^Xe[ln(ER/keVnr) — ln3] + 0.09,0} for 1 < Eji/keVnr < 3 


where Les{Eji) is the best-fit scintillation efficiency reported in Eig. 1 of Ref. [51], and ^xe is a 
nuisance parameter, first proposed in Ref. [45] to logarithmically extrapolate the scintillation 
efficiency below 3 keV. In Eq. (4.2) we also assume a Gaussian photomultiplier resolution of 
variance nupf^-p, with upmt = 0.5. Einally we adopt the detector efficiency £{Si) in Eig. 2 
of [51]. 


4.6 XENONIO 

In a second study the XENON collaboration uses the ionization signal S2 only to measure 
the nuclear recoil energy of the detected events [6] . The experimental recoil energy threshold 
for this analysis is of about 1.4 keVnr- The relation between S2 (measured in PE) and the 
observed nuclear recoil energy is 52 = Qy{Eji)Eji. We extract the function Qy{Eji) from 
Eig. 1 in [6]. The experimental exposure corresponding to these data is 12.5x1.2 kg-days. 
In this low energy threshold analysis XENONIO observed k = 23 candidate events in the 
signal region 1.4 - 10 keVnr- A large background contamination cannot be excluded [6]. We 
calculate the expected nuclear recoil energy spectrum in XENONIO as in Eq. 5.5 of [10], 
assuming a constant efficiency £{Eji) = 0.94. We also allow for = 20 background events 
of unknown origin, with ctb = 10. 


4.7 LUX 

We construct the LUX Likelihood function as for XENONIOO. The first data release of the 
LUX experiment consists of 85.3 live days of dark matter search data, corresponding to an 
exposure of 250x85.3 kg-days [7]. LUX observed 160 events with SI between 2 and 30 PE. 
Only one event is below the mean of the Gaussian fit to the nuclear recoil calibration events 
in Eig. 4 of [7]. The expected number of background events in the same region is 0.64 ± 0.16. 
We assume (Tpmt = 0.37, and the experimental efficiency in Eig. 9 of Ref. [7] divided by 2 to 
account for the 50% nuclear recoil acceptance quoted by the LUX collaboration. 

4.8 COUPP 

Using a 4.0 kg CP3I bubble chamber, the GOUPP experiment has measured k = 2, k = 3 and 
k = 8 bubble nucleations above the threshold energies 7.8 keVnr, U keVnr and 15.5 keVnr, 
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respectively [8] . We calculate the expected number of dark matter scattering events above a 
threshold energy Eth as follows [8] 







<ER,t 


(4.3) 


where e(illth) is the threshold dependent experimental exposure times the bubble detection 
efficiency, and VriER, Eti^) is the probability that an energy Er nucleates a bubble above 
£lth- The latter is given by [8] 


Vt{Er, Eth) = 1 - exp 


— Q-T 


Er - Eth 
Eth 


(4.4) 


Here we neglect scattering from Carbon (ac = 0), assume perfect efficiency for bubble 
nucleation for Iodine {a\ —>• +oo), and treat ap = ocoupp as a nuisance parameter. For the 
energy dependent exposure we assume e(7.8 keV) = 55.8 kg-days, e(ll keV) = 70 kg-days 
and e(15.5 keV) = 311.7 kg-days [8]. For the three threshold configurations the expected 
number of background events is = 0.8, hr = 0.7 and fiR = 3, respectively. 


4.9 PICASSO 

The PICASSO experiment searches for dark matter using superheated liquid droplets made 
of C4F10 [9]. In the last run PICASSO operated assuming eight different bubble nucleation 
threshold energies, namely 1.7, 2.9, 4.1, 5.8, 6.9, 16.3, 38.8 and 54.8 keV. For each energy, the 
observed count rate above threshold 77j, i = 1,..., 8, is reported in Fig. 5 of Ref. [9], including 
the associated errors (Tj. We calculate the expected scattering rate TZi in the energy range 
[i7th, -|-oo] using Eq. (4.3) with e = 1. As for COUPP, we neglect scattering from Carbon, and 
consider a-p = apicASSO as a nuisance parameter. For the eight data points TZi we assume 
the Gaussian Likelihood function: 


— In TpiCASSO = 


E 


2=1 




(4.5) 


5 Global limits and interference patterns 

We now compare the general effective theory of one-body dark matter-nucleon interactions 
to current observations, in a global statistical analysis of the direct detection experiments in 
Sec. 4. We present our results in terms of exclusion limits on the 28 coupling constants of 
the theory. Our investigation focuses on multi-interaction interference effects in the exclusion 
limit calculation. 

Multi-inter action interference effects are of two types. A first type involves pairs of dark 
matter-nucleon interaction operators. It generates terms proportional to cj'cj, with i ^ j, in 
Eq. (2.10), arising from R^, R^i, and R^-r/- Inspection of Eq. (A.l) in Appendix A 

shows that there are 7 pairs of interfering operators in the general effective theory of one-body 
dark matter-nucleon interactions, namely (OijOs), (£>4,05), (04,06), (03,09), (Oii,Oi2), 
(OiijOis), and (Oi2,Oi5). Other pairs of operators do not interfere, partially since nuclear 
ground states are eigenstates of P and CP, and partially because of their dependence. We 
refer to the remaining operators in Tab. 1 as non-interfering dark matter-nucleon interaction 
operators. 
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Table 3. Table of correlation coefficients (5.3) different from zero for the LUX experiment at = 10 
TeV. 


A second type of multi-interaction interference arises when isoscalar and isovector com¬ 
ponents of the same operator combine in Eq. (2.10) forming terms proportional to cjcj , 
T ^ t' . Any operator in Tab. 1 can produce isoscalar-isovector interference patterns of this 
type as long as off-diagonal nuclear response functions are generated in the dark matter- 
nucleus scattering. 

Let us now consider a pair of operators {Oi,Oj) and an experiment characterized by 
the signal region S. Let us also assume cj ^ 0 and cj / 0 for fixed (r, r') and {i,j) pairs, 
and all other coupling constants equal to zero. Then, contours of constant ns{Tny-,c,r]) at a 
given and rj represent ellipses in the cj — cj plane of equation 


TT T T 

an qq 


-|- 2(1, 






-I- Qjj c,- Cj = const, 


(5.1) 


where the constants oL’’ , aJJ , and ajj are given by 


= / dEn 


dTZ 


dEji 


cl=cl = 1 ; cj=cj =0 


(5.2) 


and similarly for aJJ and aJJ . 

We can therefore introduce a correlation coefficient rJJ'' to quantify the degree of interfer¬ 
ence between a pair of operators {Oi,(Dj), or between the isoscalar and isovector components 
of a given operator: 


TT _ 




(5.3) 


In Tab. 3 we list some of the rJJ' coefficients different from zero for the LUX experiment at 
rriy- = 10 TeV. 

Below we fit 77 , and c to observations globally, i.e. including all the available data, 
and simultaneously varying sets of coupling constants with rJJ' 7^ 0. 


5.1 Non-interfering operators 

We start with an analysis of the non-interfering operators. Non-interfering interaction 
operators are O7 = Sw • v"'-, Oio = ISn ■ £>13 = i{Sj^ ■ v-*-)(SAr • q/mw) and 

Oi4 = 7(8^ • q/m7v)(S7v • v-*-). In contrast to O7 and C>io, the interaction operators O13 
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Figure 1. In all panels, cyan contours represent 2D 90% confidence intervals from a fit of m^, rj and of 
the coupling constants in the legends to the direct detection experiments indicated in parenthesis. In 
each panel, colored regions correspond to the 2D profile Likelihood associated with the cyan contours. 
Yellow contours denote 2D 90% credible regions obtained by fitting and a single coupling constant 
to the LUX data. 

and Oi4 depend on both the momentum transfer operator and the transverse relative veloc¬ 
ity operator. 

All non-interfering operators contribute to the nuclear spin current through the nu¬ 
clear response operators which also appear in the theory of electroweak 

scattering from nuclei, and characterize the familiar spin-dependent dark matter-nucleon in¬ 
teraction operator O 4 . At the same time, the operator O 15 induces a nuclear spin-velocity 
current through the nuclear response operator which is specific to dark matter-nucleon 

interactions. 

The top panels in Fig. 1 show the results that we obtain fitting Cy, Cj, rj and rriy^ to 
current direct detection experiments. Results are presented in terms of 2D profile Likelihoods 
(colored regions) and associated 2D 90% confidence intervals (cyan contours) in the planes 
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Figure 2. Same as for Fig. 1 but now for the model parameters c^g, 0^4 and m^. 

— Cj (left panel) and — Cj (right panel)^. For reference, in the same panels we also 
report the 2D credible regions that we obtain fitting and a single coupling constant {cj 
in the left panel, and Cj in the right panel) to the LUX data. For single coupling constant 
fits, where marginalization is not involved, we adopt a Bayesian approach, as it requires a 
relatively small number of Likelihood evaluation. Volume and prior effects in dark matter 
direct detection are discussed in [13, 45, 52, 53]. The bottom panels in Fig. 1, and the four 
panels in Fig. 2 show the results of similar analyses which focus on the operators Oio, O 13 
and £>14. 

Comparing different panels in Figs. 1 and 2, we find that confidence intervals and profile 
Likelihoods in the rriy^ — and — cj planes are similar for a given operator Oi. This result 
is expected, in that the (0,0) and (1,1) components of the nuclear response functions W^T , 
W^T, and are comparable for a given isotope. The relative strength of the exclusion 

^More precisely, we use the decimal logarithm of the dark matter particle mass and of the coupling constants 
expressed in units of m^, with = 246.2 GeV. 
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Figure 3. In all panels, green and cyan contours represent 2D 90% confidence intervals from a 
global fit of m^, t] and of the coupling constants in the legends to all direct detection experiments in 
Sec. 4. In each panel, the colored region corresponds to the 2D proHle Likelihood associated with the 
less stringent exclusion limit at the reference value ^ 100 GeV in that figure. Yellow contours 
represent 2D 90% credible regions obtained by fitting and a single coupling constant to the LUX 
data. 


limits in Figs. 1 and 2 reflects the number of q and v;^ operators multiplying Cy, cJq, 
and c [4 in the expressions for Ig and 1^ in Eq. (2.3). In addition, it also depends on the 
relative amplitude of the nuclear response functions W^J , W^J, and integrated over 

the relevant signal regions. 

5.2 Isoscalar-isovector interference patterns 

For any interaction operator in Tab. 1, including operators that do not interfere in pairs as 
C> 7 , C>io, Oi 3 , and On, we observe isoscalar-isovector interference patterns in the rate ( 2 . 10 ). 
As an example, let us focus on the operator O^. The operator O 7 contributes to the square 
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(5.4) 


The last term in (5.4) describes a destructive interference between the isoscalar and isovector 
components of the operator £> 7 , as once integrated over a typical signal region, it gives a 
negative contribution to the total scattering rate. We observe similar interference patterns 
for all operators in Tab. 1. 

In general destructive interference effects make direct detection exclusion limits weaker. 
For instance, an otherwise excluded value of C 7 can remain compatible with observations if 
compensated by an appropriately large value of C 7 , because of the negative oc C 7 C 7 term in 
Eq. (5.4). 

For this reason, in all figures of this work single coupling constant fits to LUX data 
results in stronger exclusion limits for > 5 GeV, as in these fits destructive interference 
effects are neglected. For smaller values of m^, SuperCDMS, and CDMSlite dominate the 
exclusion limit calculation. 


5.3 Interfering operators 

Interfering operators divide into 4 independent subsets. The hrst subset consists of the 
operators Oi = I^at and O 3 = iS^v • [{^/mN) x v-*-]. The operator Oi contributes to the 
nuclear vector charge through the nuclear response operator whereas the operator 

O 3 contributes to the nuclear spin current, and to the nuclear spin-velocity current through 
respectively. The operators Oi and O 3 generate a transition probability 
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Figure 5. 2D 90% confidence intervals (cyan contours) and profile Likelihoods (colored regions) from 
a global fit of the model parameters in the legends, r) and to current direct detection experiments. 
Contours are presented in the six planes C 4 — m^, C 5 — m^, Cg — m^, C 4 — m^, Cg — m^, and Cg — 
Yellow lines represent 2D 90% credible regions obtained by htting and a single coupling constant 
to the LUX data. 
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(5.5) 


In Eq. (5.5) terms with t ^ t' describe isoscalar-isovector interference effects similar to those 
discnssed in Sec. 5.2. The term oc C 4 C 3 arises from the interference of Oi and CI 3 . Integrating 
Eq. (5.5) over a typical signal region, cancellations between different terms occur, as the 
integrated nuclear response functions WJJ and W^T , for r / r', and for t = t' , are 

negative for many isotopes. 

Because of cancellations in (5.5), numerical noise affects the exclusion limits derived 
simultaneously varying rj, m^, c^, c}, c-j, and C 3 . To circumvent this problem, while exploring 
all interference patterns, here we present exclusion limits obtained in four complementary 
ways. In a first analysis we fit m^, c^, C 3 , and r] to current dark matter direct detection 
experiments. In a second analysis, we simultaneously place limits on the constants c\ and C 3 
while fitting the same data. Our third and fourth analysis respectively consider the constants 
Cl — C 4 , and C 3 — C 3 as free parameters in the global fit. 

Eig. 3 shows the 2D 90% confidence intervals (colored contours) and profile Likelihoods 
(colored regions) resulting from the four analyses described above. To derive the contours 
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Figure 6. Same as for Fig. 5 but now for the model parameters Cg, Cg, c|, Cg and m^. 

in the panels we fit m^, r/, and the constants in the legends to observations. Results are 
presented in the — c?, — cl, — c® and — C 3 planes. For comparison, in each 

panel we also report the 2D 90% credible region that we obtain fitting rriy. and a single 
coupling constant to the LUX data. 

In the top panels of Fig 3 (corresponding to the m-^ — c? and — c} planes), we 
observe significant differences between distinct contours. As expected, limits obtained from 
single coupling constant fits to LUX data are generically stronger in the > 5 GeV mass 
region. Global fits of the c? and c\ coupling constants to current observations (cyan contours) 
instead produce the weakest exclusion limits, as in Eq. (5.5) cancellations between the term 
proportional to c[c[, with r = 0 or r = 1, and the one proportional to can be large. 
In contrast, distinct contours in the bottom panels of Fig 3 do not significantly differ, as in 
Eq. (5.5) cancellations between the term oc C 3 C 3 , with r = 0 or r = 1, and the terms oc cjcj 
or oc C 3 C 3 , with t ^ t' , are less important. 

The second subset of interfering operators consists of On = iS^ ■ (q/mAr), O 12 = 
• (Sat X V-*-), and O15 = — [S^ • (q/mAr)][(SAr x v-*-) • (q/mAr)]- The operator On contributes 
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Figure 7. Selected 2D profile Likelihoods from a fit of m^, r] and of the parameters in the legends to 
current observations. We focus on pairs of coupling constants with large correlation coefficients (see 
Tab. 3). Results are presented in the planes ch — 0 ^ 2 , 0^2 — c^g, c\i — 0 ( 2 , and 0(2 — cjg. In each panel, 
blue lines represent ellipses drawn using Eq. (5.1) with the right hand side hxed at a given reference 
value for = 10 TeV, and assuming LUX as experimental setup. 


to the nuclear vector charge through the nuclear response operator similarly to Oi. 

The operators O 12 and O 15 induce a nuclear spin current through and and a 

nuclear spin-velocity current through 

The coupling constants CJ 2 , and cjg contribute to the scattering amplitude in a way 
similar to c[ and C 3 in Eq. (5.5). Also in this case, we therefore expect cancellations in the 
scattering rate due to destructive interference between operator pairs, or between isoscalar 
and isovector components of the same operator. To study the impact of multi-interaction 
interference effects in the exclusion limit calculation, we proceed as for the operators Oi and 
O 3 . We therefore divide the coupling constants c\i, c[ 2 , and c[^ into 5 subsets: (c?^, 032 , c^g), 
(c{i, c} 2 , C 35 ) and (c?,cl), with i = 11,12,15. We then fit m^, t), and a single subset of 
coupling constants at the time to observations. Finally, from each fit we derive exclusion 
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Figure 8. Same as for Fig. 7 but now for the model parameters Cii, c\i, C 12 , and c }2 


limits on the coupling constants in analysis. This procedure allows to separate interference 
effects due to operator pairs, from those related to isoscalar and isovector components of the 
same operator. 

Fig. 4 shows the 2D 90% confidence intervals and profile Likelihoods that we find in the 
five analyses described above. Results are presented in 6 vs coupling constant planes. 
To derive the contours in the panels we fit r], and the constants in the legends to 
observations. For comparison, in each panel we also report the 2D 90% credible regions of 
a single coupling constant fit to the LUX data. In the my. — 0^2 > “ '^15 ”^x “ '^15 

planes distinct exclusion limits are comparable. In the other panels, contours differ by up 
to 1 order of magnitude in the coupling constants. As for Oi and O 3 , such differences arise 
from multi-interaction interference effects. 

There are two additional subsets of interfering operators. One consists of the operators 
O4 = • Sat, O5 = • [(q/mAr) x v-^], and Oe = [S^ • (q/mAr)][SAr • (q/mAr)]. The other 

one includes Og = • v-*- and Og = iS^ ■ [Sat x (q/mAr)]- Eq. (2.3) shows the contribution 

of these operators to the nuclear charges and currents. The operators O 5 and Og are the 
only operators contributing to the nuclear convection current (through Figs. 5 and 

6 show the 2D 90% confidence intervals and profile Likelihoods that we find analyzing the 
operators O4, O5, Oq, Og, and Og. In these figures we do not observe strong destructive 
interference patterns. 

We conclude this section with an analysis of selected pairs of coupling constants. We 
focus on largely correlated coupling constants (see Tab. 3). Figs. 7 and 8 show the 2D 
profile Likelihoods that we obtain fitting m^, r} and the parameters in the legends to current 
observations. Results are presented in the planes c^i — C 12 , C 12 — 0 ^ 5 , c\i — c\ 2 , 0^2 “ 

C]*! — and 0^2 — Ci 2 - each panel, blue lines represent ellipses constructed from Eq. (5.1) 
with the right hand side fixed at a given reference value, and assuming LUX as experimental 
setup. As expected, the red elliptical regions in Figs. 7 and 8 are in general aligned with 
the blue ellipses. However, the 2D profile Likelihood in the 0^2 — C 15 plane constitutes an 
exception. The mismatch observed in this figure is due to the fact that the blue ellipse in 
the top right panel of Fig. 7 neglects a large — 0^2 correlation, which is instead taken into 
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account deriving the corresponding red elliptical region. 


6 Conclusions 

We compared the general effective theory of one-body dark matter-nucleon interactions to 
current direct detection experiments in a global multidimensional statistical analysis. In 
this study, we included data from a variegated sample of direct detection experiments in a 
single Likelihood function. In the fits we simultaneously varied subsets of correlated model 
parameters, covering a large volume of the parameter space. We presented our results in 
terms of 90% confidence level exclusion limits on the 28 isoscalar and isovector coupling 
constants of the theory. 

We found that current direct detection experiments are able to probe all isoscalar and 
isovector coupling constants, including those measuring the strength of dark matter-nucleon 
interactions commonly neglected in this context. For instance, the interaction operators 
Oil = • (q/ruTv) and O 12 = • (Sat x v"*-), and the familiar spin-dependent interaction 

operator O 4 = • Sat are equally constrained by current data, though the former are by far 

less explored. 

We characterized the interference patterns arising in dark matter direct detection from 
multi-interaction effects involving pairs of dark matter-nucleon interaction operators, or 
isoscalar and isovector components of the same operator. Most importantly, we quanti¬ 
fied the impact of multi-interaction interference effects on the calculation of direct detection 
exclusion limits. 

We found that destructive interference effects weaken direct detection exclusion limits 
derived neglecting the superposition of different interaction operators by up to one order of 
magnitude. Destructive interference effects largely affect exclusion limits on the strength 
of the familiar spin-independent interaction Oi, and are less important for the standard 
spin-dependent operator O 4 . 
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A Dark matter response functions 


Below, we list the dark matter response functions that appear in Eq. (2.10). The notation is 
the same used in the body of the paper. In the figures, for definitiveness we assume = 1/2. 
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